Este proceso es relativamente sencillo para la base de UFOS, ya que al momento de descargar los datos de internet solo te quedas con un solo archivo el cual con crear la tabla y un copy la cargas rápido a SQL.
create table ufos (
fecha varchar,
fecha2 varchar,
ciudad varchar,
estado varchar,
figura varchar,
duracion varchar,
descripcion varchar,
reportado varchar);
\copy ufos from '/home/itam/proyectos/ufos/UFO_OK.tsv' with delimiter E'|' NULL AS 'na';
alter table ufos alter column fecha type date using to_date(fecha,'MM-DD-YY');
###Primer avistamiento por estado
select min(fecha), estado from ufos group by estado order by estado;
min | estado
------------+--------
1970-01-01 |
1970-04-15 | AB
1972-08-30 | AK
1970-03-30 | AL
1970-10-01 | AR
1970-04-15 | AZ
1970-03-15 | BC
1970-01-15 | CA
1970-01-24 | CO
1970-06-15 | CT
1992-04-15 | Ca
1974-06-06 | DC
1972-05-01 | DE
1970-01-16 | FL
1996-04-28 | Fl
1970-02-20 | GA
1970-06-06 | HI
1970-06-01 | IA
###Primer avistamiento por figura
select min(fecha), estado from ufos group by figura order by figura;
min | figura
------------+------------
1970-02-15 |
1970-05-05 | Changing
1970-09-15 | Chevron
1970-01-24 | Cigar
1970-02-20 | Circle
1970-06-30 | Cone
1997-03-22 | Crescent
1971-08-11 | Cross
1970-01-15 | Cylinder
1974-07-18 | Delta
1970-09-15 | Diamond
1970-01-01 | Disk
1996-03-15 | Dome
1970-01-20 | Egg
1970-07-22 | Fireball
###Promedio de avistamientos por mes y año
select avg(cuenta), mes
from(select count(fecha) cuenta,extract(year from fecha) anio,
extract(month from fecha) mes, estado
from ufos
group by extract(year from fecha),
extract (month from fecha), estado) as base
group by (mes)
order by (mes);
avg | mes
--------------------+-----
5.6103476151980598 | 1
4.9247121346324181 | 2
5.4217067108533554 | 3
5.1141304347826087 | 4
4.8512635379061372 | 5
4.6379944802207912 | 6
6.1336161187698834 | 7
5.9702346880366342 | 8
6.0897689768976898 | 9
5.9026143790849673 | 10
6.0802422407267222 | 11
5.6908333333333333 | 12
select avg(cuenta) promedio, anio
from(select count(fecha) cuenta,extract(year from fecha) anio,
extract(month from fecha) mes, estado
from ufos
group by extract(year from fecha),
extract (month from fecha), estado) as base
group by (anio)
order by (anio);
avg | anio
------------------------+------
1.4871794871794872 | 1970
1.3921568627450980 | 1971
1.5781250000000000 | 1972
1.5808383233532934 | 1973
1.6524064171122995 | 1974
1.7766990291262136 | 1975
1.7771739130434783 | 1976
1.7231638418079096 | 1977
1.9693877551020408 | 1978
1.7654320987654321 | 1979
1.6250000000000000 | 1980
1.5426356589147287 | 1981
1.5620437956204380 | 1982
select avg(cuenta) promedio, estado
from(select count(fecha) cuenta,extract(year from fecha) anio,
extract(month from fecha) mes, estado
from ufos
group by extract(year from fecha),
extract (month from fecha), estado) as base
group by (estado)
order by (promedio) desc
limit 10;
promedio | estado
---------------------+--------
18.9259259259259259 | CA
13.4861591695501730 |
11.6247191011235955 | FL
11.5852534562211982 | WA
9.2005649717514124 | AZ
8.9065040650406504 | TX
7.7979797979797980 | NY
7.6750000000000000 | IL
7.3729603729603730 | PA
7.0578313253012048 | OH
(10 rows)
###Estado con mayor varianza
select stddev(cuenta) desviacion, estado
from(select count(fecha) cuenta,extract(year from fecha) anio,
extract(month from fecha) mes, estado
from ufos
group by extract(year from fecha),
extract (month from fecha), estado) as base
group by (estado)
order by (desviacion) desc
limit 10;
desviacion | estado
---------------------+--------
| Ca
| VI
23.2855556245333992 | CA
15.7433042229631767 |
14.6629034962789133 | FL
13.0119278750568268 | WA
9.6515354652537361 | TX
9.4335494187814312 | IL
9.3434033517730056 | NY
9.1962414026243307 | PA
(10 rows)
Par realizar los siguientes análisis, exporte un query en csv para leerlo en R.
copy (select count(fecha) cuenta,fecha, estado
from ufos
group by fecha,estado)
to '/home/itam/proyectos/ufos/linea_tiempo.csv' delimiter ',' csv header;
Una vez filtrados los datos en Postgres los importo para analziarlos en R (la carpeta de docker está conectada con la mía de big-data).
library(dplyr)
library(ggvis)
datos<-read.csv("/home/jared/big-data/alumnos/jared275/data/ufos/linea_tiempo.csv")
##Corrijo un poco la fecha
datos$años<-as.numeric(substr(datos$fecha,1,4))
datos$años[datos$años>2015]<-datos$años[datos$años>2015]-100
datos$mes<-substr(datos$fecha,6,7)
datos$trim<-"trim_1"
datos$trim[datos$mes %in% c("04","05","06")]<-"trim_2"
datos$trim[datos$mes %in% c("07","08","09")]<-"trim_3"
datos$trim[datos$mes %in% c("10","11","12")]<-"trim_4"
Dividimos los datos en trimestres para investigar algun patrón en el tiempo de los avistamientos y notamos que el tercer trimestre es casi siempre el momento en donde más ovnis se ven, por otra parte el trimestre 1 es el periodo del año en donde menos pasan.
datos%>%
group_by(años,mes)%>%
summarise(conteo=n())%>%
arrange(años,mes)%>%
ungroup()%>%
mutate(fecha=paste(años,mes,"01",sep="-"),fecha=as.Date(fecha,"%Y-%m-%d"))%>%
ggvis(~fecha,~conteo)%>%
layer_lines()
datos%>%
group_by(años,trim)%>%
summarise(conteo=n())%>%
arrange(años,trim)%>%
ungroup()%>%
mutate(fecha=paste(años,substr(trim,6,6),"01",sep="-"),fecha=as.Date(fecha,"%Y-%m-%d"))%>%
ggvis(~fecha,~conteo)%>%
layer_points(fill=~trim)%>%
layer_lines()
Para matchear el estado con su nombre bajo una tabla de wikipedia.
library(rvest)
wiki<-html("http://en.wikipedia.org/wiki/List_of_U.S._state_abbreviations")
tabla<-wiki%>%
html_nodes(xpath='//*[@id="mw-content-text"]/table[1]')%>%
html_table()
tabla<-tabla[[1]]
tabla<-tabla[1:53,c(1,4)]
names(tabla)[2]<-"estado"
tabla_2<-left_join(tabla,datos)
tabla_2<-tabla_2[!is.na(tabla_2$fecha),]
tabla_2$Region<-tolower(tabla_2$Region)
head(tabla_2)
## Region estado cuenta fecha años mes trim
## 3 alabama AL 1 2001-03-01 2001 03 trim_1
## 4 alabama AL 1 2007-04-05 2007 04 trim_2
## 5 alabama AL 2 2011-10-21 2011 10 trim_4
## 6 alabama AL 1 2013-11-05 2013 11 trim_4
## 7 alabama AL 1 1996-06-18 1996 06 trim_2
## 8 alabama AL 1 2012-12-21 2012 12 trim_4
California concentra gran parte de los avistamientos pero si la omitimos vemos que Texas, Florida , Pensilvania y Wahington también concentran gran parte de los avistamientos a lo largo de la historia. Sin embargo, no parece que exista una correlación espacial en los avistamientos, ninguno de estos se encuentran cerca.
por_estado<-tabla_2%>%
group_by(Region,años)%>%
summarise(conteo=n())%>%
mutate(porcentaje=conteo/sum(conteo), region=Region)
map_data = ggplot2::map_data("state")
por_estado_mapa<-left_join(map_data,por_estado)
por_estado_mapa$años<-as.integer(por_estado_mapa$años)
etiquetas<-group_by(map_data,region)%>%
filter(rank(region,ties.method="first")==8)
por_estado_mapa %>%
group_by(region)%>%
mutate(conteo_tot=mean(conteo))%>%
select(long,lat,group,order,region,conteo_tot)%>%
unique()%>%
ungroup()%>%
group_by(group) %>%
ggvis(x = ~long, y = ~lat) %>%
layer_paths(fill = ~desc(conteo_tot)) %>%
layer_text(data=etiquetas,x = ~long, y = ~lat, text:=~region, fill:="white",
fontSize:=12)%>%
hide_legend(c("fill","stroke"))
por_estado_mapa %>%
filter(region!="california")%>%
group_by(region)%>%
mutate(conteo_tot=mean(conteo))%>%
select(long,lat,group,order,region,conteo_tot)%>%
unique()%>%
ungroup()%>%
group_by(group) %>%
ggvis(x = ~long, y = ~lat) %>%
layer_paths(fill = ~desc(conteo_tot)) %>%
layer_text(data=etiquetas,x = ~long, y = ~lat, text:=~region, fill:="white",
fontSize:=12)%>%
hide_legend(c("fill","stroke"))
Ahora para ver de manera un poco rudimentaria los avistamientos por estado, y ubicar si en algun momento hubo una aglomeración de avistamientos, hacemos mapas para los últimos 10 años pero vemos que se conservan las mismas proporciones.
mapas<-lapply(1:10, function(i){
por_estado_mapa %>%
filter(años==2005+i)%>%
group_by(group) %>%
ggvis(x = ~long, y = ~lat) %>%
layer_paths(fill = ~desc(conteo)) %>%
layer_text(data=etiquetas,x = ~long, y = ~lat, text:=~region, fill:="white",
fontSize:=12)%>%
hide_legend(c("fill"))%>%
add_axis("x", title =paste("año",2005+i))
})
mapas[[1]]
mapas[[2]]
mapas[[3]]
mapas[[4]]
mapas[[5]]
mapas[[6]]
mapas[[7]]
mapas[[8]]
mapas[[9]]
mapas[[10]]
La base de gdelt tiene mas sentido tenerla en una base de datos como SQL, ya que es suficientemente grande como para no poderla cargar dentro de R. El proceso de carga a la base de datos es también un poco más complicado que la de UFO. Se creó la tabla con cada uno de los campos, los nombres se obtuvieron de la págnia, y solo se cargaron los datos hasta marzo de 2013, que es el histórico que se encuentra homologado, a partir de esa fecha hay un campo nuevo en la base y para los fines de este trabajo nos quedamos hasta 2013.
Una vez creados los datos corremos un script que de manera iterativa va cargando la base.
for gdelt_file in *.csv
do
psql -d gdelt -c "COPY gdelt FROM '/home/jared/big-data/alumnos/jared275/data/$gdelt_file' WITH DELIMITER E'\t' NULL AS 'NA';"
done
Dado que nos agrada más R para hacer análisis, hemos decidido que utilizaremos las herramientas de dplyr para conectarnos con nuestra base de datos, hacer los querys que queramos e importar la información para analizarla desde R.
Algunas buenas referencias que hablan de esto las encontramos en este link y en este.
Debemos mencionar que para que este proceso fuera más fácil se instalo la base de datos en mi sistema operativo, no en el de docker.
library(dplyr)
library(RPostgreSQL)
tabla_postgres <- tbl(src_postgres(dbname="gdelt"),
,from="gdelt")
colnames(tabla_postgres)
## [1] "globaleventid" "sqldate"
## [3] "monthyear" "year"
## [5] "fractiondate" "actor1code"
## [7] "actor1name" "actor1countrycode"
## [9] "actor1knowngroupcode" "actor1ethniccode"
## [11] "actor1religion1code" "actor1religion2code"
## [13] "actor1type1code" "actor1type2code"
## [15] "actor1type3code" "actor2code"
## [17] "actor2name" "actor2countrycode"
## [19] "actor2knowngroupcode" "actor2ethniccode"
## [21] "actor2religion1code" "actor2religion2code"
## [23] "actor2type1code" "actor2type2code"
## [25] "actor2type3code" "isrootevent"
## [27] "eventcode" "eventbasecode"
## [29] "eventrootcode" "quadclass"
## [31] "goldsteinscale" "nummentions"
## [33] "numsources" "numarticles"
## [35] "avgtone" "actor1geo_type"
## [37] "actor1geo_fullname" "actor1geo_countrycode"
## [39] "actor1geo_adm1code" "actor1geo_lat"
## [41] "actor1geo_long" "actor1geo_featureid"
## [43] "actor2geo_type" "actor2geo_fullname"
## [45] "actor2geo_countrycode" "actor2geo_adm1code"
## [47] "actor2geo_lat" "actor2geo_long"
## [49] "actor2geo_featureid" "actiongeo_type"
## [51] "actiongeo_fullname" "actiongeo_countrycode"
## [53] "actiongeo_adm1code" "actiongeo_lat"
## [55] "actiongeo_long" "actiongeo_featureid"
## [57] "dateadded"
Para darnos una idea de el número de filas que tiene la base, hacemos lo siguiente dentro de la base de datos:
SELECT reltuples::bigint AS estimate
FROM pg_class
WHERE oid = 'gdelt'::regclass;
estimate
-----------
212908960
(1 row)
Como solo vamos a utilizar los datos de México hacemos un filtro para estos y nos quedamos solo con las variables que nos interesan. Incluso hacer el query en dplyr nos puede traducir el código para hacerlo en PSQL.
noticias_mex<-tabla_postgres%>%
select(sqldate,monthyear,actor1name,actor1code,actor1countrycode,actor1type1code,
actor1geo_fullname,actor2name,actor2code,actor2countrycode,
actor2type1code,eventcode,actor2geo_fullname,avgtone,actiongeo_lat,
actiongeo_long,actiongeo_countrycode)%>%
filter(actiongeo_countrycode=="MX")
show_query(noticias_mex)
## <SQL>
## SELECT "sqldate" AS "sqldate", "monthyear" AS "monthyear", "actor1name" AS "actor1name", "actor1code" AS "actor1code", "actor1countrycode" AS "actor1countrycode", "actor1type1code" AS "actor1type1code", "actor1geo_fullname" AS "actor1geo_fullname", "actor2name" AS "actor2name", "actor2code" AS "actor2code", "actor2countrycode" AS "actor2countrycode", "actor2type1code" AS "actor2type1code", "eventcode" AS "eventcode", "actor2geo_fullname" AS "actor2geo_fullname", "avgtone" AS "avgtone", "actiongeo_lat" AS "actiongeo_lat", "actiongeo_long" AS "actiongeo_long", "actiongeo_countrycode" AS "actiongeo_countrycode"
## FROM "gdelt"
## WHERE "actiongeo_countrycode" = 'MX'
Finalmente recolectamos el query y lo guardamos como un RDS.
base_bn<-collect(noticias_mex)
saveRDS(base_bn,file="/home/jared/big-data/alumnos/jared275/noticias_mex.rds")
noticias_mx<-readRDS("/home/jared/big-data/alumnos/jared275/noticias_mex.rds")
show_query(tabla_postgres%>%
filter(actiongeo_countrycode=="MX")%>%
group_by(eventcode)%>%
summarise(conteo=n())%>%
arrange(desc(conteo))%>%
mutate(conteo_per=round(conteo/sum(conteo),4)))
## <SQL>
## SELECT "eventcode", "conteo", "conteo_per"
## FROM (SELECT "eventcode", "conteo", ROUND("conteo" / sum("conteo") OVER (ROWS BETWEEN UNBOUNDED PRECEDING AND UNBOUNDED FOLLOWING), 4.0) AS "conteo_per"
## FROM (SELECT "eventcode", count(*) AS "conteo"
## FROM "gdelt"
## WHERE "actiongeo_countrycode" = 'MX'
## GROUP BY "eventcode") AS "_W1"
## ORDER BY "conteo" DESC) AS "_W2"
eventos_freq<-noticias_mx%>%
group_by(eventcode)%>%
summarise(conteo=n())%>%
arrange(desc(conteo))%>%
mutate(conteo_per=round(conteo/sum(conteo),4))
head(eventos_freq,5)
## Source: local data frame [5 x 3]
##
## eventcode conteo conteo_per
## 1 010 119395 0.0851
## 2 042 115417 0.0823
## 3 043 108903 0.0776
## 4 190 90427 0.0644
## 5 040 74095 0.0528
De los 5 eventos más frecuentes el más interesante es el 190, que se refiere a al uso convencional de fuerza por grupos organizados. Estos pueden ser policías, militares, entre otros.
noticias_mx$año<-as.integer(substr(noticias_mx$monthyear,1,4))
noticias_mx%>%
filter(eventcode=="190")%>%
group_by(año)%>%
summarise(conteo=n())%>%
ggvis(~año,~conteo)%>%
layer_bars()
Tratamos de normalizar esta gráfica considerando que es probable que se tenga más información conforme mayor sea el numero del año, por lo que dividimos el conteo de cada año sobre el número de noticias total que se colectaron en el año.
library(tidyr)
noticias_mx%>%
group_by(año)%>%
mutate(conteo_tot=n())%>%
group_by(año,eventcode,conteo_tot)%>%
summarise(conteo=n())%>%
filter(eventcode=="190")%>%
mutate(conteo_prop=conteo/conteo_tot)%>%
ggvis(~año,~conteo_prop)%>%
layer_bars()
Tenemos 2 picos temporales, uno para el año 1994, probablemente se deba al levantamiento del ejercito zapatista y en el 2010, que no sabemos que pueda ser. Vamos a analizar en que estados ocurrieron estos eventos.
library(maptools)
library(rgdal)
library(ggplot2)
mex_shp <- readOGR("/home/jared/Dropbox/Maestría CD/Estadística Computacional/mapas/estados_ligero" , "Mex_Edos")
## OGR data source with driver: ESRI Shapefile
## Source: "/home/jared/Dropbox/Maestría CD/Estadística Computacional/mapas/estados_ligero", layer: "Mex_Edos"
## with 32 features
## It has 1 fields
mex_shp@data$id = mex_shp@data$NOM_ENT
edo_df <- fortify(mex_shp, region = "id")
eventos_190<-noticias_mx%>%
filter(eventcode=="190")%>%
select(actiongeo_long,actiongeo_lat)%>%
mutate(actiongeo_long=as.numeric(actiongeo_long),
actiongeo_lat=as.numeric(actiongeo_lat))
eventos_190<-as.data.frame(eventos_190)
p <- SpatialPointsDataFrame(eventos_190, data.frame(id=1:nrow(eventos_190)),
proj4string=CRS(proj4string(mex_shp)))
proj4string(mex_shp)
## [1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
proj4string(p)
## [1] "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
res <- over(p, mex_shp)
data.frame(table(res$NOM_ENT))%>%
mutate(id=Var1)%>%
right_join(edo_df)%>%
group_by(group) %>%
ggvis(x = ~long, y = ~lat) %>%
layer_paths(fill = ~desc(Freq))
En el mapa anterior se puede ver como los estados que concentran el mayor número de noticias de este tipo son Baja California Norte, Chihuahua, Chiapas y preponderantemente San Luis Potosí.
Ahora, hay que ver como se han comportado a lo largo de 1980 a 2012.
library(RColorBrewer)
mapas_mex<-lapply(1:33, function(i){
eventos_190<-noticias_mx%>%
filter(eventcode=="190" & año==1979+i)%>%
select(actiongeo_long,actiongeo_lat)%>%
mutate(actiongeo_long=as.numeric(actiongeo_long),
actiongeo_lat=as.numeric(actiongeo_lat))
eventos_190<-as.data.frame(eventos_190)
p <- SpatialPointsDataFrame(eventos_190, data.frame(id=1:nrow(eventos_190)),
proj4string=CRS(proj4string(mex_shp)))
res <- over(p, mex_shp)
datos_hist<-data.frame(prop.table(table(res$NOM_ENT)))
names(datos_hist)<-c("Estado",paste("año",1979+i))
datos_hist
})
df_heatmap<-mapas_mex[[1]]
for(i in 2:33){
df_heatmap<-inner_join(df_heatmap,mapas_mex[[i]])
}
heat_map<-as.matrix(df_heatmap[,-1])
rownames(heat_map)<-df_heatmap[,1]
heatmap(heat_map,scale="column",Rowv=NA,Colv=NA,col=brewer.pal(9, "Blues"))
Después de revisar el histórico, el hecho de que en casi todos los años San Luis Potosí salga como el estado que más eventos de este tipo tiene, nos hace sospechar que aquellos eventos en los que no se pueda dar una ubicación más específica, ponen a San Luis como opción default, quizá porque está más o menos al centro del territorio.
También el mapa de calor nos permite ver que Chiapas y Baja California han tenido en la historia más noticias de este tipo, por otro lado Aguascalientes no ha sobresalido en ningun momento sobre este tipo de noticias, como tampoco Nayarit, Colima y Campeche.